### This script will extract all spatial data for all sites for BRT analysis

# Obtain command line arguments from .sh file

args=(commandArgs(TRUE))

# Evaluate the arguments for use in this script

for(i in 1:length(args)) 
	
	{
		
	eval(parse(text=args[[i]]))
 
	}


# Load library SDMTools

library('SDMTools')

# Establish directories

in.dir = '/home1/99/jc152199/MicroclimateStatisticalDownscale/ToAnalyse/'
out.dir = '/home1/99/jc152199/MicroclimateStatisticalDownscale/dailyspatialextracts/'
tmaxdir = '/home1/99/jc152199/MicroclimateStatisticalDownscale/AWAPASCII/RawAWAP/Tmax/'
tmindir = '/home1/99/jc152199/MicroclimateStatisticalDownscale/AWAPASCII/RawAWAP/Tmin/'
rrdir = '/home1/99/jc152199/MicroclimateStatisticalDownscale/250mASCII/Realized_Radiation/'
insoldir = '/home1/99/jc152199/GRASS/ASCIIoutput/'
setwd(in.dir)

# Read in microclimate data

micro = read.csv('micro_min_max_static_xy.csv', header=T)
micro=micro[which(micro$date<20090602),] # Get the rest of the tmax, tmin, and realrad files then delete this line and run again

# Create columns to append data to

micro$realrad = NA 
micro$AWAPmax = NA
micro$AWAPmin = NA
micro$insol = NA

# List of days

DOY.list = sort(unique(substring(micro$date,5,8)))
t.DOY = sprintf('%03i',which(DOY.list==substring(dayx,5,8)))

# Read in ASCII files relating to dayx

insol.asc = read.asc(paste(insoldir,'insoltime_',t.DOY,'.asc',sep=''))
rr.asc = read.asc.gz(paste(rrdir,'realrad_',dayx,'_WTplusbuffer_LatLong_WGS1984_9secres.asc.gz',sep=''))
tmin.asc = read.asc.gz(paste(tmindir,substr(dayx,1,4),'/tmin.',dayx,dayx,'.grid.gz',sep=''))
tmax.asc = read.asc.gz(paste(tmaxdir,substr(dayx,1,4),'/tmax.',dayx,dayx,'.grid.gz',sep=''))
	
# Append spatial data to data frame

micro$insol[which(micro$date==dayx)]=extract.data(cbind(micro[which(micro$date==dayx),c(5:6)]),insol.asc)
micro$realrad[which(micro$date==dayx)]=extract.data(cbind(micro[which(micro$date==dayx),c(5:6)]),rr.asc)
micro$AWAPmax[which(micro$date==dayx)]=extract.data(cbind(micro[which(micro$date==dayx),c(5:6)]),tmax.asc)
micro$AWAPmin[which(micro$date==dayx)]=extract.data(cbind(micro[which(micro$date==dayx),c(5:6)]),tmin.asc)

# Create a temporary data frame

t.micro = micro[which(micro$date==dayx),]
	
# Write out a .csv file of data for dayx
	
write.csv(t.micro,paste(out.dir,dayx,'.csv',sep=''),row.names=F)
	


	












